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Improving the accuracy and thus broadening the applicability of electronic density functional theory 
(DFT) is crucial to many research areas, from material science, to theoretical chemistry, biophysics 
and biochemistry. In the last three years, the mathematical structure of the strong-interaction limit 
of density functional theory has been uncovered, and exact information on this limit has started 
to become available. The aim of this paper is to give a perspective on how this new piece of exact 
information can be used to treat situations that are problematic for standard Kohn-Sham DFT. One 

O 

way to use the strong-interaction limit, more relevant for solid-state physical devices, is to define a 

new framework to do practical, non-conventional, DFT calculations in which a strong-interacting 

^ reference system is used instead of the traditional non-interacting one of Kohn and Sham. Another 

(N 

(^T) way to proceed, more related to chemical applications, is to include the exact treatment of the 

(N 

strong-interaction limit into approximate exchange-correlation energy density functionals in order 
to describe difficult situations such as the breaking of the chemical bond. 
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I. INTRODUCTION 



Density functional theory (DFT), 1 in its Kohn-Sham (KS) formulation,^ has been a real 
breakthrough for electronic structure calculations, allowing to treat systems much larger 
than those accessible to wavefunction methods. KS DFT, together with its extension to 
time-dependent (TD) phenomena (TDDFT)P made possible the theoretical study of an 
incredible huge number of chemical, physical, and biological processes. 

The key idea of KS DFT is an exact mapping^ between the physical, interacting, many- 
electron system and a model system of non-interacting fermions with the same density. 
Only one term, the so called exchange-correlation (xc) energy functional (containing all 
the complicated many-body effects) needs to be approximated. Although in principle this 
functional is unique (or "universal"), a large number of approximations have been developed 
in the last twenty years, both by chemists and physicists, often targeting different systems, 
different properties, and different phenomena. In a way, the emergence of such a "functional 
zoology" simply reflects the intrinsic difficulty of building a single general approximation able 
to recognize and capture, for each given system or process, the many-body effects relevant 
for its description. 

Despite the large number of available approximate functionals and of their successful 
applications, there are still important cases in which KS DFT can fail, which is why the quest 
for better xc functionals continues to be a very active research field (for recent reviews see, 
e.g., Refs. HHS]). For example, present-day KS DFT encounters problems in the treatment of 
near-degeneracy effects (rearrangement of electrons within partially filled levels, important 
for describing bond dissociation but also equilibrium geometries, particularly for systems 
with d and / unsaturated shells), in the description of van der Waals long-range interactions 
(relevant, for example, for biomolecules and layered materials), and of localization effects 
due to strong electronic correlations (as those occurring in Mott insulators and in low-density 
nanodevices, but also occurring in bond dissociation). These problems can hamper more 
or less severely (and sometimes in an unpredictable way) a given calculation, depending 
on their relative importance with respect to other effects that are better captured by the 
available approximate functionals. 

This work primarily aims at describing a different approach to some of the unsolved 
problems of present-day DFT, focussing on the treatment of systems with strong spatial 
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correlations. The key idea is to recognize that the non-interacting Kohn-Sham reference 
system is not always the best choice. The main idea of Kohn and Sham, which can be 
summarized as "Let's solve a model system having the same density of the physical one 
and approximate the remaining missing energy with a density functional" , can be rigorously 
generalized to model systems different from the non- interacting one of Kohn and Sham P This 
freedom can be used to choose model systems that are able to capture some of the relevant 
effects (for example near-degeneracy or strong correlations), whose computational cost is 
still low, and for which it is easier to design approximate density functionals that recover 
the missing energy. For example, in recent years this strategy has been used to address the 
problems of near-degeneracy effects and van der Waals interactions by using a model system 
with a weak long-range-only interaction (and having the same density of the physical system, 
as in KS theory). The preliminary results are so far very successfulp^^ as proved by the 
growing number of research groups that are now working on the practical implementation 
of this "short-range DFT - long-range wavefunction" (srDFT-lrWF) methodJ-^H^ 

Strong correlations, however, remain a big challenge for DFT, and in many cases are 
also beyond the reach of the srDFT-lrWF method. By "strong electronic correlation" we 
mean here the study of systems in which the electron-electron interaction largely dominates 
over the kinetic energy, creating strong spatial correlations. In such cases, it may happen 
that we need very many (billions) of Slater determinants for a proper description of the 
relevant physics, with all the natural occupation numbers becoming very small. For these 
situations both the non-interacting KS system and the weak-interacting hamiltonian of the 
srDFT-lrWF method are not the best starting point: they are not able to capture the 
physics of the system under study so that trying to describe the missing energy with an 
approximate density functional is often a daunting task (or, alternatively, the srDTF-lrWF 
method becomes as expensive as solving the Schroedinger equation for the physical system). 

In order to "visualize" this concept, Fig. [T] schematically represents the difference between 
near-degeneracy effects, characterized by the presence of few more important states with re- 
spect to the occupied KS orbitals (that can be captured with a weak-interacting hamiltonian, 
like the one used in the srDFT-lrWF method), and strong correlations, where very many 
(billions) of Slater determinants are needed for a proper description of the relevant physics 
(notice that here we are not talking about getting the energy with high accuracy, but only 
about describing the right physics: once we have a model hamiltonian which is able to do 
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that, the idea is, as in KS theory, to correct the energy with a density functional). In this 
figure levels drawn with a solid line represent the occupied KS states (labeled with "KS"), 
and dotted levels the empty ones. On the left, we have a typical near-degenerate system, 
in which few empty states strongly couple to the ground state: including them would be 
enough to describe the right physics of the system, although for an accurate energy many 
more states would be needed. On the right we have a strongly correlated system in which 
billion of states are strongly coupled to the ground state. From the point of view of the ex- 
act first-order density matrix, the first case corresponds to having some natural occupation 
numbers ni close to 1/2 (if we consider natural spin orbitals with < Uj < 1), while the 
second case corresponds to having all rii <C 1. Of course this simple, schematic, picture may 
be very different if we use a spin-unrestricted formalism to define the KS system (see also 



Sec. VI A 2), instead of a restricted one, as mostly used throughout this paper. 

Prototype systems displaying near-degeneracy effects are the Be isoelectronic series 
(where the 2s and the 2p KS levels become more and more degenerate as the atomic number 
Z increases), and the H 2 molecule along its dissociation curve, where the a g and a u KS 
energies get closer and closer as the molecule is stretched. These two simple examples are 
paradigmatic of many situations occurring in the study of chemical and physical problems, 
from heavy elements to the stretching of the chemical bond in general. A simple example 
of strong electronic correlation are low density nanodevices such as quantum dots. As the 
electronic density is lowered, spatial correlations between the electrons become stronger and 
stronger, and, as shown in Refs. |2"4"]|2"o1 for a simple model consisting of two electrons in an 
harmonic potential, all the natural occupation numbers become very small, indicating the 
presence of an infinite number of important states. In real systems studied in experiments, in 
which low-density electrons are confined at the interfaces of semiconductor heterostructures, 
this phenomenon leads, for example, to intriguing patterns in the addition energy spectra,^ 
which are suggestive of strong spatial correlations and have never been fully explained. 

Of course, in general there are very many different physical situations which need a huge 
number of Slater determinants to be described, and many corresponding ansatz wavefunc- 
tions, models or methods that can do that, each one being able to capture different phys- 
ical phenomena. Typical examples are the density-matrix renormalization group (DMRG) 
method, the Laughlin wavefunction, the unrestricted Hartree-Fock plus symmetry restora- 
tion wavefunction, and dynamical mean field theory. 
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FIG. 1: Schematic illustration of the difference between near-degeneracy effects, in which few more 
important states with respect to the Kohn-Sham occupied levels are needed in order to capture 
the right physics, and strongly-correlated systems, which need billions of Slater determinants. The 
first case is usually characterized by the presence of natural occupation numbers ni close to 1/2, 
while the second case often corresponds to natural occupations that are all very small. In this 
figure levels drawn with a solid line represent the occupied KS states, and dotted levels the empty 
ones. On the left, we have a typical near-degenerate system, in which few empty states strongly 
couple to the ground state: including them would be enough to describe the right physics of the 
system, although for an accurate energy many more states would be needed. On the right we have 
a strongly correlated system in which billion of states are strongly coupled to the ground state. 

The main object of this paper is to review and discuss the perspectives of a new way to 
deal with the case of strong spatial correlations in a DFT framework. For a given iV-electron 
system with density p(r), we construct, in a mathematical rigorous way, a model system 
consisting of N electrons having the same density p(r) and maximum possible correlation 
between the N electronic positions. We call this model system the "strictly correlated 
electron" (SCE) model, and we use it as a complementary alternative to the KS ansatz 



5 



for DFT. We also propose simple approximate density functionals to recover the difference 
between the energies of the physical system and of the SCE model, following the same ideas 
used in KS DFT. The SCE model is able to capture the infinitely many Slater determinants 
needed to describe strong spatial correlations, and, as we shall see in the next sections, it 
is the natural counterpart of the KS ansatz. It also provides a rigorous lower bound for the 
exact exchange-correlation functional of KS DFT, simply because the electrons cannot be 
more correlated than the SCE state in a given one electron density p(r). 

The paper is organized as follows. After reviewing the basics of DFT in Sec. [TTJ, in order 
to emphasize the analogies and the differences between the usual KS DFT and the "SCE 



DFT", we parallel, throughout Sees. |IIIj]VI[ the two approaches. Thus, Sees. |IH|]VI| contain 
a KS part, which quickly reviews the main formalism pertinent to the KS ansatz, and a 
SCE part, which explains how the same concepts can be generalized using the SCE model 



as a reference system. In Sec. VII we report first applications of the SCE-DFT method 
to few-electron quantum dots at low density. Although, as previously mentioned, bond 
dissociation can be viewed as a near-degeneracy effect (which can be described by the weak 
interacting hamiltonian of the srDFT-lrWF method or, e.g., by density matrix functional 
theory^ZEHl or by a mixture of Hartree-Fock and Hartree-Fock-Bogoliubov methods^), it 
is also characterized by strong spatial correlations between the electrons involved in the 



stretched bond, whose physics can be captured by the SCE limit. In Section VIII , thus, 
we discuss possible ways to include the exact information contained in the SCE limit into 
functionals useful for chemical applications, with emphasis on bond dissociation. The last 



Sec. IX is devoted to conclusions and perspectives. 



II. THE HOHENBERG-KOHN FUNCTIONAL AND ITS BASIC PROPERTIES 

We begin by defining the problem and reviewing the basic properties of the Hohenberg- 
Kohn functional. 

We generally consider here systems of N interacting electrons, bound by a given external 
potential v(r) in D-dimensional space (re R D ). The corresponding Hamiltonian, 

N 

H a [v]=f + aV ee + J2 v ( r i), (2- 1 ) 

1=1 
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with the universal operators of the kinetic energy, 

h 2 N d 2 



2m f-f dr 2 

i=i 1 

and the Coulomb repulsion between the electrons, 

2 N 1 r 

^4E^ (2-3) 

has four independent parameters: the particle number N, the spatial dimension D, the 
R D — > R function v = v(r) of the external potential, and a tunable dimensionless interaction 
strength a > (which will be set to its realistic value a = 1 at the end). Unlike a and t> , 
the parameters N and D will not be indicated explicitly in our notation. 
Due to the Ritz principle, the ground-state energy of H a [v] is given by 

E a [v} = mui(^\H a [v}\^), (2.4) 

where the condition ^ — > N addresses all (normalized) spin-^ fermionic wave functions in 
D-dimensional space, 

V = -9(ri,...,r N ;(Ti,...,(T N ), (2.5) 
with Yi G R D and spin variables cij. A considerably simpler function is the particle density, 

2 



p(r) = N ^2 j d D r 2 ~.d D r N ty(r,r 2 ,...,r N ;cri,...,cr N 



(2.6) 



<T 1 ,...,cr N 



which is normalized according to J d D rp(r) = N. In terms of this function as the variable, 
the universal Hohenberg-Kohn (HK) functional of DFT is defined a d 30 * 31 * 

F a [p] = mm(q\f + aV C e\^) > (2.7) 
where the condition ^ — > p now addresses only those fermionic A^-electron wave functions 



\I> that are, via Eq. (2.6), associated with the same given particle density p = p{v). Here, 
"universal" means that F a [p] does not depend on the parameter v = v(r). [It does, however, 
depend on the spatial dimension D and on the particle number N — J d D rp(r).] If the 
functional F a [p] was known explicitly in terms of the density p, the ground-state energy of 



Eq. (2.4) could be obtained by a considerably simpler minimization procedure, 

E a [v) = min ^F a [p] + J d D rp(r)v{r)} (2.8) 



where the condition p — >■ N now addresses all (non-negative) density functions p(r) that are 
normalized to the same given particle number N. Eq. ( 2.8[ ) is called the (second part of 
the) HK theorem [the first part being the statement that the external potential v(r) in the 
Hamiltonian of Eq. ( |2.1[ ) is unambiguously fixed by its ground-state density p(r)]. 

Introducing a Lagrangian multiplier \i to account for the condition p — > N (and writing 
F a= i[p] = F[p\), we obtain from Eq. (2.8) the Euler equation 



SF[p] 
5p{r) 



+ v(r) 



(2.9) 



to be solved for the wanted density function p(r). Since F[p] is not known explicitly in terms 
of the density p, the crucial problem of DFT is to find approximate ways of treating F[p] 
and its functional derivative 5F[p]/5p(r). 

Clearly, the complexity of the many-body problem is hidden in the HK functional F a [p\. 
An equivalent functional is 



F p [p\ = min(^|/3T + V ee \V) = /3F 1/P \p]. 



(2.10) 



Since a minimizing wave function here at the same time minimizes Eq. (2.7) for the inter- 



action strength a = 1//3, the parameter /3 may be dubbed the "interaction weakness". 



For a given density p and interaction strength a in Eq. (2.7), let ^ a [p] be a minimizing 
wave function. With T a [p] = a [p\\f\^ a [p\) and Vee\p\ = a [p}\V ee \^ a [p\) we have 



F a [p]=T a [p] + aV e { a \p]. 



(2.11) 



We make here the usual assumption that ^ a [p] depends smoothly on the parameter a. (This 
assumption may break down, e.g., for a uniform electron gas at low density going through 
a ferromagnetic phase transition). Then, F a [p] is differentiable with respect to a and, due 



to the minimum property, Eq. (2.7), the Hellmann-Feynman theorem implie; 



da 



F a \p] = {* a \p]\V~\*a\p]). 



In particular, we can write Eq. (2.11), in terms of the universal functional^ 

Vi a) [p] = ^F a [p], T a [p] = F a [p] - a-^FM 



An immediate consequence of Eq. (2.12) is the coupling- const ant integra 132 ' 

F 1 \p]-F \p]= I daV^\p\. 



(2.12) 



(2.13) 



(2.14) 



In an analogous way, the corresponding formula for the functional Fp[p] is obtained, 

F 1 [p]-F [p] = fdpfplp] (2.15) 
J o 

(notice that F x \p) = F x \p)). Here, fp[p] = (V p[p]\f\V p[p)} where Vp[p) = ^ a =i/p[p] is 



a minimizing wave function in Eq. (2.10), Tp[p] = T a= yp[p]. Substituting j3 



a 1 , we 



obtain™ 3 

F 1 [p]-F [p}= r^T a \p]. (2.16) 

We define a density p to be ground-state-(a, v )-representable if there exists a single- 
particle external potential t> a [p](r) (whose existence is not always granted^) such that p is 
a ground-state density of the Hamiltonian 



N 



H a [p] = f + aV ec + ^2v a \p](r i ). (2.17) 

i=l 

In this case, ^ a [p] is a ground state of H a [p\; the corresponding ground-state energy, 



E a [p] = F a [p] + j d D rv a [p}(r)p(r), 



(2.18) 



however, can be degenerate. 
Similarly, the Hamiltonian 

N 



H p \p\ = pf + V ee + ^2v p \p](Ti), vp = /3v a=1/ p (2.19) 



i=i 



has the ground state ^p[p] = *&i/p[p] and the ground-state energy 

Ep[p) = Fp[p] + j d D rvp[p](v)p(r). (2.20) 

III. ZERO AND STRICT COULOMB CORRELATION 
A. Non-interacting electrons (NIE) 

The usual Kohn-Sham system corresponds to the non-interacting limit a = of the HK 
functional -F Q [p], 

F [p] = lim \Fp[p] = min(tt|T|*) = T s [p]. (3.21) 
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Being a ground state of the non-interacting Hamiltonian H a =o\p], the minimizing wave 
function *&o[p] — ^oo[p] = ^ mE [p] is, in most cases, a single Slater determinant of N spin- 
orbitals <f>i(r,a) which obey the Kohn-Sham (KS) single-particle Schrodinger equations 

{ " 2^T V2 + w o[p](r)}^(r,(7) = e^a). (3.22) 

Consequently, T s [p] is the kinetic energy of N = J d D rp(r) non-interacting electrons in a 
given ground-state density p = p(r). By construction, the KS potential u [p](r) is such that 
the orbitals reproduce the given density, 

^|<Mr,<x)| 2 = p(r). (3.23) 

i.a 

Implicitly, in terms of these orbitals (rather than explicitly in terms of the density p itself), 
T s [p) is given by 

T *M = ^E / ^|V«Mr,a)| 2 . (3.24) 

Non-interacting electrons (NIE) have zero Coulomb correlation. For example, N = 2 
such electrons in a given density p(r) have opposite spins and occupy the same spatial 



orbital ip(r) = y |p(r) (the situation can become more complicated if the corresponding 
KS potential has a degenerate ground state, something that rarely happens for N — 2). 
When their two positions are measured simultaneously, the results r x and r 2 are completely 
uncorrelated - when only the partial result ri is noticed while the result r 2 is ignored or 
hidden, its probability distribution is rigorously independent of the particular value of ri. 
In this case, the expectation of V ee is given by 

V2>W = e 2 / i°n f i°r 2 W£lM£j)£ = l Ulfi] (N = 2) , (3.25) 

J J | r l — r 2| * 

with the explicit density functional of the Hartree energy, 

V[A^h°r[i°r>*^l. (3.26) 



2 . 

If the electrons were repulsive bosons (b), an arbitrary number iV of them could occupy the 



same orbital ip(r). In this case, Eq. (3.25) would be generalized to 14e [p] = '[p] where 

V^[p} = ^U[p] (bosons). (3.27) 
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For N > 3, however, non-interacting electrons must occupy two or more different orbitals. 
Consequently, their positions can no longer be completely uncorrelated. This effect is some- 
times called Pauli correlation, since it is not caused by a true repulsive (Coulomb) force 
between the electrons, but merely by the Pauli principle. As a result, the true value of 



Ke°^[p] is for iV > 3 lower than the bosonic value of Eq. (3.27), 

K ( e 0) [p] = U[p]+E x [p}<^U[p}. 
The exchange energy E x [p] < is another implicit density functional, 



(3.28) 



EM 



e 

T 



/ j m\,m s I 



d D r I d u r 



D 



>r( r )^(r)^(r')^(r' 



r — r 



(3.29) 



with <pi(r, a) = ipi(r)Xm i ( cr )- m Eq. (3.28), the equal sign, implying E x \p] = — -hU[p), holds 
for N <2, while E x [p\ < -j^U[p] for N > 3. 



B. Strictly correlated electrons (SCE) 

In the case a > 0, the Coulomb repulsion between the electrons is turned on in the 



Hamiltonian H a [p] of Eq. (2.17). Now, the ground state ^ a [p] has, in addition to Pauli 
correlation (for N > 3), also true Coulomb correlation which is caused by a repulsive force 
which lowers the value of V&t [p] as a grows. Here we consider the extreme limit a — > oo of 
infinitely strong repulsionj^SES which we call the "strictly correlated electrons" (SCE) limit, 

lim -F a [p] = F [p] = min(^|Vy^> = V e f E [p}. (3.30) 

The functional Kg CE [p] is the natural counterpart of the KS non- interacting kinetic energy 
T s [p] and was first addressed about ten years ago, | 39 | 4C | but only treated in an approximated 
way, using physically motivated models.®^ Only recently V^ CE [p] and the square | [p] \ 2 = 
l^o[p]| 2 — |^I/ SCE [p]| 2 of the corresponding minimizing wave function have been treated 
exactly in Ref. H21 where the interested reader can find more mathematical details. In 
the following, we summarize the basics of the SCE solution, describing the physics that is 
captured by V^ e CE [p}. 

Ke CE H corresponds to the lowest possible value of the expectation of the electron-electron 
repulsion in a given density p(r). In other words, the functional Ke° E [p] defines a classical 
problem with a given smooth density. Thus, in contrast to ^ a [p] for finite a < oo, the 
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limiting wave function \Poo[p] does no longer depend on the spin variables <Ti, <7at, and, 
since the limit is classical (even if it is an unusual classical problem because of the constraint 
of the smooth density), we can only determine l^oofp]! 2 , which, in terms of the spatial 
variables ri,...,rjv, is no longer a regular function, but rather a Dirac-type distribution, 
describing electrons with strictly correlated positions. In practice, this means that the N 
results r-j e R D (i = 1,...,N) of a simultaneous measurement of all electronic positions in 
the distribution l^oofp]! 2 are no longer independent of each other, but strictly related via N 
so-called co-motion functions fj(r), 



fi(r) (z = l,...,iV; fi(r) 



(3.31) 



In other words, the position ri of one electron fixes the positions (i > 1) of all the others. 
The co-motion functions obey the group properties^ 

(3.32) 



{f 1 (f„,(r)),...,f^(f n (r))} = jf^r),...,^)} (n = 1, N), 



so that Eq. (3.31) does not conflict with the symmetry postulate on a wave function for 



identical fermions. Moreover, as the position of one of the electrons determines the positions 
of all the others, the probability of finding one electron at position r in the volume element 
d D r must be the same of finding the i th electron at position fj(r) in the volume element 
d D fi(r). This means that all the co-motion functions for a given iV-electron density p = p(r) 
must satisfy the differential equation^ 



p{U{r))d u m=p{r)d D r (x = 1 JV), 

whose initial conditions are fixed by making the corresponding V°F E [p], given by 

1 - 6« 



(3.33) 



K S e CE M 



e 
~2 



N 

E 



dPr P(r) 



N \fi(r) 



(3.34) 



minimum.^1. Thus, similarly to the N single-particle orbitals 4>i(r, a) in the NIE Kohn-Sham 
state, the co-motion functions f$(r) are fixed by the given density function p = p(r)^ 



Equation (3.34) should be viewed as the counterpart of Eq. (3.24) which, also implicitly, 
represents the density functional T s [p] = T NIE [p] for the non- interacting kinetic energy in 
terms of the orbitals 4>i(r, a). The latter represent the counterpart of the co-motion functions 
fj(r) in Eq. (3.34). The counterpart of Eqs. (3.28) and (3.29) for the functional <\p], in 
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FIG. 2: The probability distribution for the electron-electron distance ryi for the H~ anion calcu- 
lated with a very accurate wavefunction for the physical system, with the "exact" Kohn-Sham (KS) 
Slater determinant (built from a very accurate density), and with the strictly correlated electron 
(SCE) construction. All quantities are in Hartree atomic units. 

contrast, is the limit a — > oo of T a [p], which, as we shall see later, must be treated with some 
care since it diverges but still yields a finite "first-order" correction to the energy functional 

The two functionals T s [p] and Kg CE [p] define two different and complementary model 
systems in which the one-electron density is the same. A simple way to grasp the very 
different physics captured by the two model systems is to look at the probability density 
P( r i2) of finding two electrons at a distance r 12 . As an example, in Fig. [2] we report 
this probability P(riz) f° r t ne H~ anion calculated using a very accurate wavefunction for 
the physical system (see Refs. and references therein), using the Kohn-Sham non- 

interacting Slater determinant (constructed from the same accurate density), and using the 
SCE construction (see also Ref. H5|) . The three probabilities P(r 12 ) correspond to three 
systems having the same one-electron density, that is, the same probability to find one 
electron at r in the volume element d 3 r. As we see from Fig. the probability distribution for 
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the electron-electron distance is very different: in the KS system there is a higher probability 
of finding the two electrons close to each other than in the physical system, in which there 
is Coulomb repulsion. In the SCE state, the two electrons never get closer than a certain 
distance tq 4.2 a.u., and they avoid each other as much as possible without breaking the 
constraint of being in the given one-electron density. 



C. Density scaling 



For a given density p = p(r), we consider the usual continuous set of scaled densities 
Pa(t), 

p A (r) = A D p(Ar) (A > 0). (3.35) 
The prefactor X D guarantees that J d D r p\(r) = J d D r p(r) for all A > 0. 



As the orbitals 0i(r, a) solve the KS equations Eq. (3.22) and yield in Eq. (3.23) the 
density p(r), the scaled orbitals (f)[ X \r,a) = X D ^ 2 (f)i(Xr,a) yield the scaled density p\(r), 

N 

££k (A W)l 2 = PA(r), 



(3.36) 



i=l a 

and solve the modified KS equations 

h 2 



2m, 



V 2 + ^o[pA](r)}# ) (r, ( r) = efV! A) (r,a), 



(3.37) 



where vo[p\](r) = A 2 vo[p](Ar) and ef^ = A 2 ej. Therefore, Eq. (3.24) implies^ 

T s [p A ] = A 2 T s [p]. 

For completeness, we note 

that™ 

u\p x ] = \u[ P ], e x [ Px ] = xem k ( c 0) M = AK { c 0) [p]. 



(3.38) 



(3.39) 



Similarly, as the co-motion functions fj(r) solve the SCE equation (3.33) for the density 
p(r), the scaled co- motion functions 

1 



A 



fi(Ar) 



(3.40) 



solve the corresponding equations for the scaled density p\(r). Consequently, Eq. (3.34) 
implies the scaling behavior 

K S e CE N = AV e f E [p]. (3.41) 
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We notice that the HK functional has a more involved scaling behavior,^ 

F[p x ] = \ 2 F a [p] (a = A" 1 ), (3.42) 



which is an immediate consequence of Eq. (2.7) with Eqs. (2.2) and (2.3). Thus, for finite 



a (0 < a < oo), we could, without loss of generality, confine ourselves to the case a = 1. 



IV. WEAK AND STRONG COULOMB CORRELATION 



Dropping the subscripts a and the superscript (a) in Eq. (2.11), we now address the 
realistic situation with interaction strength a = 1, 

F[p}=T[p] + V cc [p}. (4.43) 

Here, F[p] = F a= i[p], etc. For N > 2, the two contributions on the right-hand side obey the 
relations 

T[p] = (* a=1 [p]\f |*«=i|>]) > T s [p] = min(*|f |*> > 0, (4.44) 

* >/' 

V ee [p] = (y a =M\V C e\y a =M) > V c f E [p] = min(^|Vy^> > 0. (4.45) 

(In the trivial case N — 1, of course, we have T[p] = T s [p) > and V ee [p] = V^ e CE [p] = 0.) 
These inequalities hold, since the realistic wave function \l/ a =i [p] is significally different from 
each one of the two minimizing wave functions on the right-hand side, \I/a=o[p] = \I/ NIE [p] 
and ^/oo[p] — ^ [p], respectively. While the latter ones are characterized completely by iV 
single-particle orbitals <fti(r,cr) or, respectively, by N co-motion functions fj(r), the realistic 
wave function \1/ Q= i[p] is mathematically much more involved. Describing electrons with 
finite Coulomb repulsion, it has neither zero nor strict, but rather some finite Coulomb 
correlation, a situation which is much harder to describe mathematically. 



The non-interacting kinetic energy T s [p] in Eq. (4.44) can be considered as the zero- 
point kinetic energy resulting (by the uncertainty principle) from the spatial confinement 
of non-interacting electrons in the density p = p(r). For interacting electrons (a = 1), this 
zero-point energy is increased by Coulomb correlation, since one such electron, due to the 
repulsion by the other ones, has less effective space available than a non-interacting one 
(a = 0) within the same given density p = p(r). Consequently, the resulting difference, 

T c [p]=T[ P ]-T s [p]>0, (4.46) 
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is called kinetic energy due to correlation. [We note in passing that, as a — > oo grows beyond 



its realistic value a — 1, this zero-point energy grows indefinitely, see Eq. (4.50) below. 



On the other hand, increasing Coulomb repulsion (a — > oo) lowers the expectation of 
the operator V^c (which is a measure for the average inverse distance |r — r'| _1 between two 



electrons in the state The second inequality, Eq. (4.45), expresses the fact that this 



lowering is maximum in the limit a — > oo of strict correlation, while it is lesser in realistic 
systems with a = 1 and finite correlation. Therefore, following Ref. |36l the difference 



V d [p] = V ee [p] - V e f E [p] > 



(4.47) 



is called decorrelation energy' 
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Combining the fundamental scaling law of Eq. (3.42) with the expressions in Eq. (2.13), 



one finds the individual scaling properties of the functionals T[p] and Ke[p], 

(a = \- 1 ), 



Tip,} = \ 2 T a [p], 

V ee [p X ] = \Vee\p\ 



(4.48) 



in contrast to Eqs. (3.38) and (3.41). From section III, we know the finite limits 



tonTa\p] = TU lim K ( C Q) [P] = V®\p], lim V^[p] = V e f E [p}. 



In addition, we have the divergent limitP^I 



a^oo: T a [p) ^ T ZP [p]a 1/2 + 0(a°), 



(4.49) 



(4.50) 



where Tzp[p] is the leading coefficient of the expansion describing zero-point oscillations 
of strictly correlated electrons about the SCE limit Consequently, the high-density limit 



A — > oo : 



(HDL) of Eq. ( gj§ reads 

T[p x ] A 2 T s [p] = T s [p x ], 
V ee [ Px ] XV^[p] = V^[p x } 

In the low-density limit (LDL), in contrast, we have 

T[p x \ -> A 3 / 2 T ZP [p] = T ZP [p A ] 



(4.51) 



A : 



(4.52) 



Here, we have used Eqs. (3.38), (3.39), (3.41), and the relation Tzp[pa] = A 3//2 Tzp[p] from 
Ref. El 
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Now, we see that the kinetic energy T[p\] in the HK functional 

F\p x ]=T\p x ] + V ee \p x ] (4.53) 

becomes dominant and approaches its non-interacting value T s [pa] in the HDL (A — > oo), 
while in the LDL (A — > 0), the potential energy Ke[PA] becomes dominant and approaches its 
strictly correlated limit Ke° E [PA]- Therefore, we call an electron system with given ground- 
state density p weakly correlated (WCOR), when T[p] ^> V ee [p] or, more precisely, 

F[p}ZT[p]>T s [p] »Kc[p] (4.54) 

and strongly correlated (SCOR), when V ee [p] ^> T[p] or, more precisely, 

F[p] % VM Z Kf E [p] » T[P]- (4-55) 



V. APPROXIMATING THE HK FUNCTIONAL 



A. Exchange-correlation (xc) and kinetic-decorrelation (kd) energies 



When the single-particle orbitals </>i(r, a) of Eq. (3.24) and the co-motion functions fj(r) of 



Eq. (3.34) can be constructed rigorously for any given density p = p(r), the functional T s [p] 
and K| CE [p] can be treated exactly. Consequently, there are two natural ways of partitioning 
the HK functional F[p}. The usual one of Kohn and Sham, 



F[p] = T s [p] + Eg[p], Eg[p] = T c [p] + VM, 



(5.56) 



treats T B [p] exactly, and looks for an approximation to the remaining contribution 



Since F[p) = F x \p] and T s [p] = F [p), Eq. fl2T4| ) now reads 

El[ P ]= tdaV^[p\. (5.57) 

The KS DFT scheme works well for weakly and moderately correlated systems (WCOR). 
For SCOR systems, where F[p\ is dominated by Ke CE [pL better results should be obtained 
by partitioning the HK functional as 



F[ P ] = Ke CE [p] + EM, E kd [p] = T[p] + V d [p}; 



(5.58) 
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with V^ CE [p] to be treated exactly and £?kd[p] to be approximated. Eq. (2.15) now reads 

da 



E kd [p] = / d(3T p [p\ 



Q 



T a [p}. 



(5.59) 



The natural counterpart of this so-called kinetic-decorrelation (kd) energy^ 2123 [p] is 
the xc-Hartree energy -E^.[p] of Eqs. ( 5.56|5.57 ). This functional is usually written as 



E*IP}=E XC [P] + U[p], 



with the functional of the exchange-correlation (xc) energy, 



(5.60) 



e xc [p] = K ( e 0) [p] - u[ P ) + VM - vi%) + T c [pl 



(5.61) 



E*\p] E c [p] 

where we have introduced the correlation energy E c [p\. An equivalent representation is 



Note also that 



E xc [p] = (Tip] - T s [p]) + (V e eM - U[ P ] 

E kd [p] = T s [p] + E xc [p] - (V c f E [p] - U[p] 



(5.62) 



(5.63) 



B. Local-density approximation (LDA) for E xc [p] and -Ekd[p] 

A simple approximation to the functional E xc [p] or, equivalently, = E xc [p] + U[p] 

is the local-density approximation (LDA), 



£ X L C DA [P] = / d D rp(r)4?)(r s (r)). 



(5.64) 



As a function of r, the dimensionless local density parameter r s (r) is given by 



r s (r) 



p(r)B 



D 



1/D 



P( r ) 



B D r s (v)^ 



(5.65) 



where Bp is the volume of a .D- dimensional ball with radius a# = h 2 /m e e 2 . E.g.: B 3 = ^a%, 
B 2 = ita 2 B . The crucial quantity in Eq. (5.64) is e x ^\r s ), the xc energy per particle in the 
D- dimensional uniform electron gas with (uniform) density p = {B^t^)^ 1 . 

The functions e x ^\r s ) for D = 2, 3 are not known analytically, but accurate parametriza- 
tions of numerical Quantum Monte Carlo (QMC) data are available. In the case D = 2, 
the data and parametrization of Attaccalite et alW^ are nowadays widely used. For D — 3, 
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popular parametrizations of the Ceperley and Alder QMC dataware the ones of Vosko, Wilk 
and Nusair^ and of Perdew and Wang.^D Remarkably, the function e^(r a ) can be interpo- 
lated accurately between its high- (r s <C 1) and low-density (r s 3> 1) limits, almost without 
relying on any QMC input at all.^ Finally, for D = 1 parametrized QMC data of the ground 
state energy of a uniform electron gas with regularized electron-electron interaction are also 
available pl 
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Given E^ A [p], a corresponding LDA for i?kd[p] is readily obtained from Eq. (5.63) 



E^ A [P} = j d D rp(v)e£\r s (r)), (5.66) 
with the kd energy per particle in the D-dimensional uniform electron gas, 

^\r s )=ti D \r s ) + ^(r s )- a -^. (5.67) 

' s 

The non-interacting kinetic energy ts(r s ) per particle in the uniform electron gas (in units 
of 1 Ha = e 2 /as = m e e 4 /h 2 ) is known analytically, 

*W-ifi±O^L^-!+£ (, 68) 

and the coefficient aj^ determines the Madelung energy (in units of 1 Ha), 

ajj = -1.1061, ajj = -0.89593. (5.70) 

a (D) 

The Madelung energy exactly corresponds to the thermodynamic limit (number of 
particles and volume going to infinity with the particle density kept fixed) of V^p E [p]/N in 
a uniform electron gas (with the usual cancellation between the Hartree term, the electron- 
background and the background-background interaction energies). Thus, as in KS theory, 
the LDA is uniquely defined as the approximation that makes the method exact in the limit 
of uniform density. 

C. Exact first-order approximation for E xc [p] and i?kd[p] 

In KS DFT the exact first-order approximation for E xc [p] is the exchange energy of 



Eq. (3.29), which, as said, is an implicit functional of the density through the KS orbitals. 
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The "first-order" approximation for £>kd[p] corresponds to zero point (ZP) oscillations 
around the SCE minimum.^ The proof that this is indeed the exact first-order correction is 
rather lengthy and the interested reader can find all the details in Ref. H3 

Basically, in the SCE limit the total potential energy of a classical configuration 

e 2 

E pot (r 1 ,...,v N ) = I — - — | + ^ v sce[p](^i) , (5.71) 

Kj J i 

where t>scE[p](r) is the external potential associated with the density p at zero kinetic en- 
ergy, is constant on the D-dimensional subspace Q = {fi(r), . . . , fjv(r)} of the full ND- 
dimensional configuration space^ and is expected to have a minimum with respect to varia- 
tions perpendicular to Qq, implying that its Hessian has D eigenvectors with null eigenvalue 
and ND — D eigenvectors with positive eigenvalue 0J^(r) at every point on Qq. 47 In terms 



of these eigenvalues, the small /3 and the large a expansion of Tp\p] defined after Eq. (2.15) 



and T a [p) of Eq. (2.11) read 



hm7>[p] = r 1/2 Tz P [p]+O(/3 ) 

/3-s>0 

lim T Q [p] = a 1/2 T ZP [p] + 0(a°), 



(5.72) 
(5.73) 



with 



Tzp[p\ 



d D r ^ N J2^1 



N 



(5.74) 



Thus, as anticipated in Sec. IV, in the strict correlation limit the kinetic energy grows 



indefinitely. However, both Eqs. (5.72) and (5.73) when inserted in Eq. (5.59) yield the 
finite result 

„ , . ND-D ( v 

Eg\A = 2T ZP [p] = dPr # V (5.75) 



N 



At=l 



which is the SCE counterpart of the exact exchange energy of Eq. (3.29) for KS theory. The 
energy ^^[p] has a highly non trivial functional dependence on p, so that its functional 
derivative is not easily accessible. 
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VI. EXACT TREATMENT OF T B [p] OR V^ e CE [p] 



A. The Kohn-Sham approach (exact T s [p]) 



1. Spin-restricted formalism 



With Eq. (5.56) for the HK functional F[p], the Euler equation Eq. (2.9) reads 



8p(r) 

with the electrostatic potential 



$[p](r) + v xc [p](r) +v(r) = fi, 



$[p](r) 



5p(r) 



r — r 



of the density p(r) and the xc potential, 



SE xc [p] 



(6.76) 



(6.77) 



(6.78) 



Sp(r) 

When the approximation E^[p] used for E xc [p] is an explicit density functional, the cor- 
responding functional derivative f^[p](r) = 5E^[p]/5p(r) can be evaluated for any given 
density function p(r). 



By varying the density through variations of the orbitals, Eq. (6.76) for interacting elec- 



trons is formally equivalent to the corresponding equation for a system of non-interacting 
electrons in the KS effective external potential 



VKs[p)(r) = $[p](r) + 

Thus, the KS orbitals satisfy the equations 

h 2 



r) + v(r) 



2m f 



V 2 + «KsM(r)}&(r,<7) = ef s <Mr,a), 



(6.79) 



(6.80) 



which have to be solved self-consistently with Eq. (3.23). 



Since the exchange-correlation functional must be approximated in practice, one obtains 
an approximate ground-state energy for the physical interacting system, E^ p = T s [p] + 
(E^[p} + U[p}) + Jd D rp(r)v(r). 



Employing in Eq. (6.80) the exact quantum-mechanical operator of the kinetic energy, the 



functional T s [p] is treated exactly here. Consequently, this approach works well in the case 
of WCOR systems when T s [p] is the dominant contribution to F[p\. For SCOR systems, in 
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contrast, we will analyze in the next section a complementary approach based on the exact 
treatment of V^ e [p]. Before doing so, however, we briefly review the widely used spin-DFT 
(or unrestricted Kohn-Sham) formalism. 



2. Spin-unrestricted formalism 

In practical calculations the spin-DFT version^ of KS DFT is widely used. Although 
the Hoehenberg-Kohn functional only depends on the total density p(r), in spin DFT one 
introduces the functional T s [p^,p±\, 

T s [p t , Pi \ = min (*|T|*>, (6.81) 

which corresponds to the kinetic energy of a non-interacting system having given spin den- 
sities Pt( r ) an d P4-( r )) with 

2 

*(r,r 2 , ...,r N ;a, a 2 ,...,a N ) , (6.82) 



p CT (r) = N J d D r 2 ...d D r N 

(T2,...,<7jv 



and p^ + pi = p. The functional T s [p^,pi\ can be used to decompose the HK functional as 
F[p]=T s [p h p i ] + U[p]+E xc [p hPi } + J d D rv(r)p(r), (6.83) 



where E xc [p^, pj is defined as the correction needed to make Eq. (6.83) exact. The idea is 
to have a non-interacting system with the same spin densities of the true, interacting, one. 
This constraint defines two effective potentials t>Ks,t[p]( r ) an d ^Ks,4.[p]( r ), and two sets of 
orbitals such that £V |</> ij0 -(r)| 2 = p a (r). 

Notice that we have (for the exact functionals evaluated at the exact density and spin 
densities) T s [p t ,pJ > T s [p], £ xc [p t ,pJ < E xc [p], and T s [p t ,pJ + E xc [p t ,pJ = T s [p] + E xc [p\. 
Using the spin-unrestricted KS reference system instead of the restricted one allows to mimic 
some correlation effects, similarly to the spin-unrestricted Hartree Fock method. 



B. The SCE approach (exact V^ e CE [p\) 

The non-interacting functionals T s [p] and T s [p^, pj require a self-consistent procedure for 
their calculation. This is because the density (or the spin densities) is determined by the 
KS orbitals by the simple equation ^ |0j(r)| 2 = p(r), while determining the orbitals from 
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the density requires a highly non-trivial procedure (for which very many different numerical 
techniques have been proposed in the last years, e.g.J^HSU). 

The construction of the complementary functional V^ CE [p] for strictly correlated electrons 
for a given density p(r) can be simpler, because the density determines the co-motion func- 



tions fj(r) via the differential equations (3.33). In other words, in the SCE case it is easier 



to determine the co-motion functions from the density than to determine the density from 
the co-motion functions. In particular, Vg| CE [p] has been directly constructed for spherically 
symmetric densities, 42 while algorithms to solve the SCE equations in the general case are 
under study: a very promising way to proceed is to exploit the similarity between the SCE 
problem and mass transportation theory.^ 

The problem of calculating V^i CE [p] can be reformulated as^21 

K S e CE [p]=min / |^( ri ,r 2 ,...,r w )| 2 ^^^, (6.84) 

V>~W ~ \Ti - Tj\ 

where \ip\ 2 is the spatial part of the many-electron wavefunction. As said, in fact, in the 
SCE case, the electrons are strongly distinguished by their relative positions, so that the spin 
state (or more generally, the statistics) does not play a role The functional V^ CE [p] is thus 
the same as the spin unrestricted functional V^ e CE [p^, pj (of course with p-j- + p± = p) . This 
means that also the exact kinetic and decorrelation functional i?kd[p] is the same in the spin 
restricted and spin unrestricted formalism. However, when we deal with approximations for 



£/kd[p] this might not be true. In Sec. VII B we will compare the results for a quantum dot 
with three electrons obtained by using the local spin density functional E^ B [p^, pj with 
those from the LDA functional. 

Since the co-motion functions can be constructed from the density, in the SCE ap- 
proach we can obtain the many-electron energy by directly minimizing the expression 



F[p] + f d D rv(r)p(r) with respect to the density function p(r), according to Eq. (2.8). To this 



end, the HK functional F[p] must be partitioned as in Eq. (5.58) where an approximation 
i^kdH is required for the functional -Ekd[p], 

E ap [v] = mm {v c f E [p} + Eg\p] + J d D rv{v) p{v)) . (6.85) 

Unlike the KS equations, this approach should be particularly suitable for SCOR systems for 
which the HK functional is dominated by Kg CE [p]. In such cases, the density is dominated 
by strong spatial correlations rather than by the quantum mechanical shells. In practical 
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calculations, the minimization of Eq. (6.85) can be carried out by expanding the density on 
a suitable basis set or by using a grid. A simple example of such a calculation is reported 
in the next Sec. [VII Al 

Another equation that the minimizing density must satisfy can be obtained by varying 
the energy with respect to p(r): 



5E[v\ 5V*x\p\ | SEM | 



(6.86) 



5p(r) Sp(r) <5p( r ) 
where p is the chemical potential. Although the functional V^ CE [p] depends on the density 



in a rather complicated way via the co-motion functions [see Eq. (3.34)], its functional 
derivative t>scE[p]( r ) = — ^fe^ satisfies the classical equilibrium equation^ 



N 



Vw S CE[p](r) = 



i=2 



(6.87) 



which has a very simple physical meaning: the potential t>scE[p](r) must compensate the net 
force acting on the electron in r, resulting from the repulsion of the other N — 1 electrons 
at positions fj(r). The one-body potential t>scE[p]( r ) is the counterpart of the KS effective 



potential of Eq. (6.79) and corresponds to the Lagrange multiplier for the constraint \I/ — > p 



in the minimization of Eq. (3.30 ). Thus, another possibility to solve the SCE-DFT equations 



is to look for the density p(r) that satisfies Eqs. (6.86), (6.87) and (3.34). This last way 
to proceed, however, raises some questions about the uniqueness of the solution, questions 
that will be addressed in future work. 



VII. SCE-DFT APPLIED TO FEW-ELECTRON QUANTUM DOTS 

In this Section we report preliminary applications of the SCE-DFT method on simple 
quantum dots models with few electrons. 

Quantum dots are nanodevices in which the motion of electrons is quantized in all 
three dimensions through the lateral confinement of a high-mobility modulation-doped two- 
dimensional electron gas in a semiconductor heterostructure (for a review, see, e.g.j^). Be- 
cause the confinement of electrons in these "artificial atoms" can be varied at will, they 
have become a playground in which the basic physics of interacting electrons can be largely 
explored and theoretical models can be tested. The number of confined electrons can vary 
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from a few to several hundred, with smaller numbers of electrons becoming increasingly 
technologically important in nandevices such as the single-electron transistor. 

In quantum dots the correlation effects between electrons need to be considered carefully 
because the external confinement can become much weaker than in real atoms, where the 
independent electron model with mean-field theories usually gives good results. As the con- 
finement strength is lowered, the mutual Coulomb interaction becomes gradually dominant. 
The physics of this regime can be thus much better captured by SCE-DFT than by tradi- 
tional KS-DFT. Indeed, KS DFT has proved useful for studying quantum dots in the weakly 
correlated regime (e.g.J^^j), while the medium and strongly-correlated regime, and in par- 
ticular the cross-over from the Fermi liquid behavior to the Wigner-crystal-like state, has 
only been accessible to wavefunction methods, e.g., configuration interactio n 1 59 1 64 1 6 ^ (only 
for very small dots), Quantum Monte Carlo (e.g.j^SS) or unrestricted Hartree-Fock plus 
symmetry restoration.^! Here we explore with SCE-DFT the regime of weak confinement 
(strong correlation), where state-of-the-art KS-DFT breaks down. 

We thus consider a simple quantum-dot model consisting of N electrons in two dimensions 
(2D) laterally confined by a parabolic potential: 

f-2 N 2 N N 2 N 

i=l i=l j=i+l J i=l 

where m* is the effective mass and e the dielectric constant. 

For now we only analyze single dots for which we obtain circularly symmetric densities, 
p(r) = p(r). In this case, the problem of determining V^ CE [p] can be separated into an 
angular part and a radial part/ 42 - The distance r from the center of the dot of one of the 
electrons can be freely chosen, and it then determines the distances from the center of all the 
other N — 1 electrons via radial co-motion functions fi(r), as well as all the relative angles 
Oij(r) between the electrons.^ The radial co-motion functions fi(r) can be constructed as 
follows.^ Define an integer index k running for odd iV from 1 to (N — l)/2, and for even N 
from 1 to (JV - 2) /2. Then 



/2*(r) 



N~\2k- N e (r)) r < a 2k 
N- 1 (N e (r)-2k) r>a 2k 



/2*+i(r) = <^ 1 (7.89) 
N- 1 (2N-2k-N e (r)) r > a N . 2k , 
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where a,- = N„ 



N e {r) 



27r xp(x) dx, 



(7.90) 



and N e 1 (y) is the inverse function of N e (r). For odd N, these equations give all the needed 
N — 1 radial co-motion functions, while for even N we have to add the last function, 



f N (r) = N; 1 (N-N e (r)). 



(7.91) 



The relative angles 0^ (r) between the electrons can be found by minimizing numerically the 
electron-electron repulsion energy 2^2i>j[fi( r ) 2 + fj(r) 2 — 2fi(r)fj(r) cos fly] -1 ' 2 . The radial 



co-motion functions of Eqs. (7.89)-(7.91 ) satisfy Eq. (3.33) for 2D circularly symmetric p, 



2tt fi(r)p(fi(r)) \f-(r)\dr = 2nrp(r) dr } 



(7.92) 



and, together with the minimizing angles %(r), yield the minimum expectation of V et 



Physically, the solution of Eqs. ( 7.89 )-( 7.91 ) makes the N electrons always be in N different 



circular shells, each of which contains, on average in the quantum mechanical problem (at 
a = 1), one electron. In the SCE limit, the electrons become strictly correlated, and all 
fluctuations are suppressed (see, e.g.,^): the space is divided into N regions, each of which 
always contains exactly one electron. 



A. The case N = 2 



In this case the minimizing angle is always Ouir) = it and there is only one co-motion 
function given by 

f 2 (r) = N-\2-N e (r)), (7.93) 

with /2(/2( r )) — T i thus ensuring the equivalence of the two electrons. 

We switch to effective Hartee units ih — 1, a* B — ^ctB = 1, e = 1, m* = 1), and we 
define f(r) = .MO, so th & t 

Vim ^ r^rf^ = f ^ r -M- y ,7,4) 

where we have used the fact that, since the electrons are indistinguishable, integrating from 
to oo is equivalent to integrate N times from to a\ = A* r ~ 1 (l). This is a characteristic of 
the SCE limit: the space is divided in N equivalent regions, so that to calculate the energy 
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we only need to treat one of them. In a way, the SCE limit seems to become more "local", 
a characteristic which may prove very useful if we deal with approximations. However, we 
also have to keep in mind that, although for an exact evaluation of V^ CE [p] we need indeed 
only one of the N equivalent regions, in order to find how to divide the space in those N 
regions we need often to perform a classical minimization over the whole space. This will 
become clearer in the next example with N = 3 electrons. 

The exact "first-order" or zero-point energy is, in this case, given by 

EII\P\ = I 1 dr7rrp(r)Mr)+w 2 (r)], (7.95) 

with 



o 



= \i f^+Tf^ (7 - 96) 

rf(r) (r + f(r)) 



W2(r) = t/ + (7 . 97) 

1 '-/'(r) (r + /(r)) 3 ^ ' 

In Ref. |37] we have evaluated the energy functional _E ap [t>] = V^ CE [p] + + 
J d D rv(r)p(r) using the exact input densities from Ref. [7TJ and we have compared the results 
with standard KS-LDA ones (notice that for two-dimensional electronic structure calcula- 
tions LDA is still the most widely used functional). At this postfunctional level we have found 
that, as expected, for large values of the confining parameter u (corresponding to higher 
densities) the KS LDA result is superior to the SCE-DFT. However, as u becomes smaller 
(which corresponds to lowering the density and thus approaching the strongly-correlated 
regime), the SCE-DFT results with its approximations for -E^dM become better and better, 
highly outperforming KS-LDA. These results are summarized in Fig. [3j where we report 
the absolute % error on the total energy as a function of the confinement parameter oj for 
KS-LDA and for SCE-DFT with Eg\p] = (curve labeled SCE), with E^[p] = E^ A [p] of 



Eq. gggfr (SCE-LDA), and with E^[p] = E^[p] of Eq. ( [5/75] (SCE-ZP). For the ground 
state energy of the 2D electron gas (which defines the LDA functional) we have used the 
data and parametrization of Attaccalite et al. 4 ^ We see from Fig. [3] that for uj < 0.007 the 
SCE-ZP result is the most accurate. The much simpler SCE-LDA is also very reasonable in 
this regime, reducing the error of KS-LDA by a factor 5-10. 

The next step is to perform self-consistent SCE-DFT calculations, in which the density 
is determined by minimizing the energy functional. Here we report very preliminary results 
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FIG. 3: The absolute % error on the total energy as a function of the confinement parameter cj 
made by the functional E^[v] = V£P E [p] + E^[p] + / d D rv{r)p(r) with Egftp] = (SCE), with 



^kd A W of Eq- ( 5 - 66 > (SCE-LDA), and with = ^[p] of Eq. (|5?75|) (SCE-ZP). 



The results obtained with standard KS-LDA are also reported. In this figure all calculations are 
done at the postfunctional level only. 



obtained by parametrizing the density with a set of N g gaussians: 

p [p] {r) = C- 1 (^c l e-^ , (7.98) 

where {p} denotes the set of the 2N g variational parameters i = l,...N g }. The 

constant C ensures that p{r) is normalized to AT = 2 electrons, and the functional form 
guarantees that p(r) > everywhere. As an example, here we consider two cases with small 
confining parameter, u = 0.0072846 and u = 0.00221088, for which we find that N g = 3 
gaussians are enough to accurately reproduce the exact density (when the fitted densities 
are inserted in i? ap [f] the error with respect to the energy obtained with the exact densities 
is ~ 0.01%). We consider only the simple SCE-LDA functional and perform the direct 
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FIG. 4: Radial densities for N = 2 electrons in a two-dimensional model quantum dot for two 
different values of the confining parameter uj. The exact values^ are compared with the results 



obtained by the direct minimization of the energy functional SCE-LDA of Eq. (7.99). Effective 
Hartree atomic units are used. The corresponding total energies have relative errors, respectively, 
of 5.4% and 4.4%. 



minimization 

E ap [v] = mm {V e S e CE [p W ] + C[P W ] + / d D rv(r)p {p} (r)} (7.99) 

with respect to the parameters {p}. This way of proceeding is probably not the best one both 
in terms of efficiency and accuracy, but the aim here is only to show a proof of principle. 
Better procedures are currently under study. The minimizing densities are compared in 
Fig. [4] with the exact ones obtained from the solution given in Ref. [7H Although the 
densities obtained are quite reasonable, it is evident that the LDA approximation for the 
functional £kd[p] has a tendency to give densities that are too diffuse. The total energies 
obtained in this way are quite accurate, with errors of 5.4% (for u = 0.0072846) and 4.4% 
(for u = 0.00221088), corresponding, respectively, to absolute errors of 3 mH* and 1 mH*. 



B. The case N = 3 



In this case we have two co-motion functions, /2(f) and /a(r), and two relative angles 
that have to be minimized numerically for each value of the distance r e [0, cq] of one of 
the electrons from the center of the dot. Notice that if, say, electron 1 is in the circular 
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shell < r < a\, then electron 2 is in the shell a\ < fair) < and electron 3 is in 
a 2 — fs{ r ) < 00. Thus, even if we only need to compute the minimizing angles for r e [0, ai], 
we explore the whole space where p(r) 7^ through the positions of the other N— 1 electrons. 

The quantum dot with N = 3 electrons is also a useful example to discuss the spin state 
in the framework of SCE-DFT. Accurate wavefunction methods, in fact, (see, e.g. 

,1313) find 

that the ground state for the N = 3 dot with oj < 0.05 is fully spin polarized. As discussed in 



Sec. VI B the functional V e g CE [p], being essentially classic, is independent of the spin state. 
The exact functional £7kd[/°] should thus be the same as the exact functional E kd [p^, pj, when 
the exact density or the exact spin densities are used. When constructing approximations, 
however, one could obtain better results with E kd [p^, pj, as in KS-DFT. 

Here we consider only the SCE-LDA and SCE-LSD functionals, and we apply them at the 
postfunctional level using as input the Diffusion Monte Carlo densities from Refs. 16611681 We 
study the values u = 0.01562, 0.005 and 0.001, which already lie in the regime where KS-LDA 
orbitals become difficult to obtain (notice that the KS-LDA results of Fig. [3] were obtained 
at the postfunctional level, using the exact densities as input). As said, we explore the two 
options -E'k < P A [p] and ££f D [p-|-, pj for which we use the parametrization of the 2D electron 
gas energy of Attaccalite et alW^ This functional is based on accurate Diffusion Monte Carlo 
(DMC) data predicting a weakly first order transition from the unpolarized gas to the 
fully polarized state at r s « 26. Even if the existence of this transition has been recently 
questioned in Ref. [721 we stick here to the original Attaccalite et al. parametrization. Since 
the densities involved are quite low, corresponding often to r s > 26, the correct definition 
(within the chosen parametrization) of the LDA functional consists in taking in each point 
of space the ground state energy of the electron gas with the same density, i.e., 

£kd DA H = J d D r p(r) {e kd (r s (r), ( = 0) 9 (25.56 - r s (r)) + e kd (r.(r), ( = 1)9 (r s (r) - 25.56)} 

(7.100) 

where r s (r) = (7rp(r))~ 1//2 , C = (Pt — Pl)/P> an d is the Heaviside step function. For the 
values of the confinement parameter u considered here (for which the ground state of the 
dot is fully polarized), instead, the "exact" LSD functional (i.e., the one which has not only 
the exact local density in each point of space, but also the exact local spin densities) is 



^kTM = J d D rp(r)e kd (r s (r)X = 1). (7.101) 
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TABLE I: Relative % errors on the total energy of a model two-dimensional quantum dot consisting 
of 3 electrons confined in an harmonic potential u ex t( r ) 



7jW 2 r 2 . Columns as follows: SCE are 



the results obtained by setting £?kd[p] = 0, SCE-LDA are those obtained by using E^ A [p] of 



Eq. (7.100), and SCE-LSD are those obtained by using E^ D [p^,0] of Eq. (7.101). 



SCE SCE-LDA SCE-LSD 



0.01562 -15.1 3.4 3.9 
0.005 -10.6 3.6 3.7 
0.001 -6.7 2.8 2.8 



In Table [T] we report the % errors on the total energies (with respect to the DMC energies) 
obtained with the two functionals. We also show the results corresponding to -E^dH = 0, 
labeled "SCE". We see that the quality of the two local approximations is rather good, with 
the LSD results slightly worse than the LDA ones for u = 0.01562 and u> = 0.005. This 
is due to the fact that for these values of the confining parameter u, r s (r) is often still 
smaller that 26, so that a lower energy is obtained by considering the true ground state of 
the electron gas. At w = 0.001, we have r s (r) always greater than 26 so that LDA and 
LSD become the same. In other words, the SCE-LDA functional predicts a transition to 
the fully polarized state at a much lower uj with respect to the one predicted by accurate 
wavefunction methods. This transition in the SCE-LDA method entirely depends on the 
delicate physics of the 2D uniform electron gas, and it is thus questionable in view of the 
latest results of Ref. [721 

This simple example shows that the next step for the construction of functionals useful 
for SCE-DFT is probably by considering simple exchange models, which would allow to 
distinguish between different spin states, generalizing to nonuniform densities what has 
been done for the uniform electron gas in Ref. [731 

VIII. IS THE SCE LIMIT RELEVANT FOR CHEMICAL APPLICATIONS? 

The results of the previous Section suggest that the SCE formalism can have an impact 
on solid-state devices involving electron gas in low dimensional systems (quantum wires, 
dots, point contacts, etc.), in the low-density, strongly-interacting regime, where traditional 
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KS DFT is not of much use. It is however less evident whether the SCE limit could be also 
relevant for applications in chemistry. 

If we consider the simplest chemical system, the H 2 molecule, we see that, as we stretch 
the chemical bond, the energy and physics of the system is exactly described by the SCE 
limit, as electrons in a stretched bond have strong spatial correlations (see, e.g., Fig. 11 of 
Ref. |74J). This feature is very interesting and promising, since the stretching of the chemical 
bond is one of the typical situations in which restricted KS-DFT encounters problems, being 
unable to describe the strong correlation occurring between the electrons involved in a single 
or in a multiple bond. Thus, the SCE limit contains useful exact information that is usually 
missed by state-of-the-art (restricted) KS-DFT. However, when we deal with real chemical 
systems the situation is different from that of the simple H2 molecule, since only the electrons 
involved in the stretched bonds are strongly correlated. The SCE limit applied to the whole 
system would give much too low energies, producing serious over correlation. In other words, 
we cannot expect the SCE-DFT scheme to work for chemistry, where often both the orbital 
description and strong spatial correlation are important at the same time. 

What we could do, instead, is trying to include the exact information contained in the SCE 
limit into approximate exchange-correlation functionals. Attempts in this direction have 
been done in the past, leading to the construction of the interaction-strength- interpolation 
(ISI) functional! 39 * 75 ! As shown in Eq. (5.57), the exchange-correlation energy of KS-DFT is 



given by (in this section we use Hartree atomic units) 

E xc [p]= [ daVW[p]-U[p]. (8.102) 
Jo 

Since the functional V^[p] approaches the SCE limit as a — > 00, the idea of the ISI 
functional is to construct the a-dependence of W a [p] = Vje [p] — U[p] by interpolating 
between the a — >• (exchange energy and second-order Gorling-Levy perturbation energy^ 

£ C GL2 H), 

W a ^ ip] = V c ^r%] - U[p] = E x [p] +2«£ c GL2 [p] + 0(« 2 ), (8.103) 
and the a — > 00 limits (SCE plus ZP oscillation^ 



W^ip] = Vjr^lp] - U[P\ = V e s e CE [p] - U[p\ + %M + 0(a"«) q > - (8.104) 

Ja 4 



However, this way of proceeding leads to serious size-consistency errors. The size-consistency 
problem of the ISI functional is related to the fact that the interpolation is done on the global 
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quantity W a [p]. Moreover, when the ISI was first proposed an exact treatment of the SCE 
limit was not available, so that the functional relied on physical approximations for the SCE 
and ZP energies] 41 * 75 * 

As a possible way out, the exact solution of the SCE limit, now available, makes accessible 
not only global, but also local quantities. This new access to local quantities could be used to 
construct local interpolations along the DFT adiabatic connection, restoring size consistency 
(for critical reviews on the size-consistency issue in DFT see alsd 77 * 78 *) . We thus rewrite 



Eq. (8.102) in terms of an energy density w a (r; [p]), 

E xc [p] = j d D rp{v) J daw a (r;[p}), (8.105) 



with 



J d D rp(r)w a (r; [p]) = W a [p] = v£>\p\ - U[p}. (8.106) 
The idea is then to use the energy densities w a (r; \p\) in the a — > and a — > oo limits, 



describing locally the quantities of Eqs. (8.103 )-(8. 104), in order to construct an interpolation 



for the a— dependence of w a (r; [p]). Since the energy density w a (r; [p]) is not uniquely 
defined, we must use the same gauge for the weak and and the strong-interaction limits. A 
very reasonable and physical choice would be the gauge defined by the exchange-correlation 
hole, 

M*M) = \ J d D u^^, (8.107) 

where u = r 2 — ri, u = |u| and the exchange-correlation hole p" c (r,u) is simply related to 
the pair density P^(ri,r2) obtained from the wavefunction 



P 2 a (r u r 2 ) = N(N-l) J d D r 3 . . .d D r N \^ a (r 1 ,a 1 , . . .r N ,a N ) 



I 2 , 



(8.108) 



<T1,...,<T2 



p« c (r, u) = — J £ (P 2 a (r, r + u) - p(v)p(v + u)) . (8.109) 

The a — )■ limit of w Q (r; [p]) is thus the exchange energy density defined in the gauge of the 
exchange hole, for which one could use the exact exchange hole or a good approximation, 
e.g., the one of Becke and Roussel.^ The a — > oo limit of w a (r; [p]) is exactly given by the 
SCE solution, which is already defined in the gauge of the exchange-correlation hole (see 
also Ref.H5l), 



•Wr, H) = 1 t |f . ( ; ) :^ (r)| - / ^(r + u). (8.110) 

l,j=l 



33 



Much more difficult is to have a local expression for the next leading terms, both for a — > 



and a — > oo, defined in the same gauge. The zero-point term of Eq. (5.75), which determines 
how the a — > oo limit is approached to orders a~ 1//2 , is, in fact, expressed in a gauge which 
is not the one of the exchange-correlation (xc) hole. The Gorling-Levy perturbation theory 
is also difficult to define locally in terms of the xc-hole gauge. 

Routes to define and calculate the local next leading terms will be pursued in future 
work. For the ZP term, one could actually directly calculate the pair-density associated to 
the 0(a' 1 ^ 2 ) wavef unction, S3 and produce the exact exchange-correlation hole in this limit. 
For the a — y leading correction, one should probably use different correlation-strength 
indicators than the GL perturbation theory. A very promising route could be the one 
described by Becke in Ref. EDI which considers the local normalization of the exact exchange 
hole as an indicator of strong non-dynamical correlation. 

The main message of this Section is that the SCE limit contains useful exact information 
for critical situations in Chemistry such as stretched bonds. However, one has to be able 
to use this exact information locally, where it is needed. This direction of research will be 
pursued in future work. 



IX. CONCLUDING REMARKS 



The strong-interaction limit of density functional theory, exactly solved in the last three 
years, contains useful physical and chemical information, typically missed by standard Kohn- 
Sham DFT. In this paper we have outlined some paths to fully exploit this piece of exact 
information, with the aim of broadening the applicability of DFT for electronic structure 
calculations in solid-state physical devices and in chemical systems, addressing fundamental 
issues of standard KS DFT. 

The mathematical structure of the strong-interaction limit of DFT has been uncovered in 
Refs. I4*2ll45|l471 However, solving the relevant equations for a general density in an efficient 
way is still an open problem, which will be addressed in future work, exploiting the formal 
similarity with mass transportation theory.^ 

Another line of research for future work is based on the fact that the strictly correlated 
problem defined by the strong-interaction limit of DFT provides a physical, rigorous, lower 
bound for the exact exchange-correlation functional of standard Kohn-Sham DFT, a feature 
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which may be exploited for the construction of approximate functionals.^ 

The calculation and study of energy densities in the strong-interaction limit of DFT will 
also provide useful information to be included into approximate functionals, and will be the 
object of future work. 
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